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Abstract. In the pattern matching approach to imaging science, the 
process of "metamorphosis" is template matching with dynamical tem- 
plates [28] . Here, we recast the metamorphosis equations of [28] into the 
Euler-Poincare variational framework of [14] and show that the meta- 
morphosis equations contain the equations for a perfect complex fluid 
[13]. This result connects the ideas underlying the process of metamor- 
phosis in image matching to the physical concept of order parameter 
in the theory of complex fluids. After developing the general theory, 
we reinterpret various examples, including point set, image and density 
metamorphosis. We finally discuss the issue of matching measures with 
metamorphosis, for which we provide existence theorems for the initial 
and boundary value problems. 



1. Overview 

Pattern matching is an important component of imaging science, with 
privileged applications in computerized anatomical analysis of medical im- 
ages (computational anatomy) [2, 4, 12, 23]. When comparing images, the 
purpose is to find, based on the conservation of photometric cues, an optimal 
nonrigid alignment between the images. In this context, diffeomorphic pat- 
tern matching methods have been developed, based on this principle, and on 
the additional goal of defining a (Riemannian) metric structure on spaces of 
deformable objects [7, 24]. They have found multiple applications in medical 
imaging [16, 3, 11, 30, 10, 29], where the objects of interest include images, 
landmarks, measures (modeling unlabeled point sets) and currents (model- 
ing curves and surfaces). These methods address the registration problem 
by solving a variational problem of the form 

(1) Find CLTguiin ( d(id : g)^ 1 + Errorterm(^.77.f em p, ntarg) ) 
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over all diffeomorphisms g, where n temp and nt arg are the compared objects 
(usually referred to as the template and the target), (g,n) i-> g.n is the ac- 
tion of diffeomorphisms on the objects and d is a right-invariant Riemannian 
distance on diffeomorphisms. This problem therefore directly relates to geo- 
desic in groups of diffeomorphisms, namely to the EPDiff equation [1, 18, 14], 
and the conserved initial momentum that specifies the solution has been 
used in statistical studies in order to provide anatomical characterizations 
of mental disorders [22, 31]. 

One of the issues in problems formulated as (1) is that the error term 
breaks the metric aspects inherited from the distance d on diffeomorphisms. 
This model has an inherent template vs. target asymmetry, which is not 
always justified. With the purpose of designing a fully metric approach 
to the template matching problem, Metamorphoses have been introduced 
in [20], and formalized and studied in [28, 27]. They provide interesting 
pattern matching alternatives to (1), in a completely metric framework. 
In this paper, we pursue the following twofold goal: (i) Provide a generic 
Lagrangian formulation for metamorphoses, that includes the Riemannian 
formalism introduced in [28], and (ii) Study a new form of metamorphoses, 
adapted to the deformation of measures. 

We start with (i), for which our point of view will be fairly abstract. 
We consider a manifold, N which is acted upon by a Lie group G: N 
contains what we can refer to as "deformable objects" and G is the group 
of deformations, which is the group of diffeomorphisms in our applications. 
We will review several examples for the space N later on in this paper. 

Definition 1. A metamorphosis [28] is a pair of curves (g t , rft) £ G x N 
parameterized by time t, with go = id. Its image is the curve n t G N defined 
by the action nt = gt-Vt- The quantities gt and rjt are called the deforma- 
tion part of the metamorphosis, and its template part, respectively. When 
rjt is constant, the metamorphosis is said to be a pure deformation. In 
the general case, the image is a combination of a deformation and template 
variation. 

In [28] , metamorphoses were used to modify an original Riemannian met- 
ric on N by including a deformation component in the geodesic evolution. 
In this paper, we generalize the approach to a generic Lagrangian formu- 
lation, and apply the Euler-Poincare variational framework [14] to derive 
evolution equations. More specific statements on these equations (for exam- 
ple, regarding the existence and uniqueness of solutions) require additional 
assumptions on G and the space N of deformed objects. In the second 
part of this paper, we will review the case in which N is a space of linear 
forms on a Hilbert space of smooth functions, which will allow us to define 
metamorphoses between measures. 

The next section provides notation and definitions related to the general 
problem of metamorphoses. 
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2. Notation and Lagrangian formulation 

We will use either letters rj or n to denote elements of N, the former being 
associated to the template part of a metamorphosis, and the latter to its 
image. 

The variational problem we shall study optimizes over metamorphoses 
(gt,r]t) by minimizing, for some Lagrangian L, 

(2) / L(g t ,g t ,7] t ,f]t)dt 

Jo 

with fixed boundary conditions for the initial and final images no and n\ 
(with nt = gtVt) an d go = idc (so only the images are constrained at the 
end-points, with the additional normalization go = id). 

Let denote the Lie algebra of G. We will consider Lagrangians defined 
on TG x TN, that satisfy the following invariance conditions: there exists 
a function I defined on q x TN such that 

L(g,U g ,r),£ v ) =Z{U a g~ 1 ,gr),gt,r ) ). 
In other terms, L is invariant by the right action of G on G x N defined by 
(g,v)h = {gh,^ 1 ^). 

For a metamorphosis (g t ,rj t ), we therefore have, letting ut = gtg^ 1 , fH = 
g t r]t and v t = g t rj t , 

L{gt,gt,Vt,Vt) = i(u t ,n t ,u t ). 

The Lie derivative with respect to a vector field X will be denoted Cx- 
The Lie algebra of G is identified with the set of right invariant vector fields 
U g = ug, u G TidG = g, g £ G, and we will use the notation C u = Ljj. 

The Lie bracket [u, v] on q is defined by 

and the associated adjoint operator is &d u v = [u,v]. Letting I g (h) = ghg~ l 
and Ad v g = £ v I g (id), we also have ad u v = £ u (Ad„)(id). When G is a group 
of diffeomorphisms, this yields ad u v = duv — dv u. 

The pairing between a linear form I and a vector u will be denoted (l | uj . 
Duality with respect to this pairing will be denoted with a * exponent. 

When G acts on a manifold N, the o operator is defined on TN* x N and 
takes values in q*. It is defined by 

[5 on| u) = — (<5| un). 

3. Euler Equations 

We compute the Euler equations associated with the minimization of 

l 

l{ut,n t ,v t )dt 

with fixed boundary conditions n$ and n\. We therefore consider variations 
5u and uj = Sn. The variation bv can be obtained from n = grj and v = gr) 



L 
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yielding h = u + un and Co = Su + mo + Sun. Here and in the following of 
this paper, we assume that computations are performed in a local chart on 
TN with respect to which we take partial derivatives. 
We therefore have 

rl / .SB , , §£ 



. \ (St \ (SI 



Co t - u t io t ~ Su t n t ) dt = 0. 



)) 



10 \^Su 
The Su term yields the equation 

S£ SI 

— + — o nt = 0. 
on dz^ 

(Note the abuse of notation StjSv G T(TN)* is considered as a linear form 
on TiV by (o~£/oV| z) := (<5^/<5i^| (0,2;)).) For the to term, we get, after an 
integration by parts 

d_st_ 
dis~v~ + ' 

where we have used the notation 
(3) 



,5£, 
Su ' 



Si 


SI 


Ju ~ 


Sn 


[u * 


SI | 



We therefore obtain the system of equations 



(4) 



St St 

+ o n< = 



3 5£ 



<9t oV U * oV o~n 



nt = u t + w t nt 



Note that + |^ o n is the momentum arising from Noether's theorem 
for the considered invariance of the Lagrangian. The special form of the 
boundary conditions (fixed no and n\) ensures that this momentum is zero. 



4. EULER-POINCARE REDUCTION 

An equivalent system can be obtained via an Euler-Poincare reduction 
[14]. In this setting, we make the variation in the group element and in the 
template instead of the velocity and the image. We let & = Sgtg^ 1 and zu t = 
gtSrjt- From this, we obtain the expressions of Su, Sn and Su. We first have 
Sut = £t + [£t,«t]; this comes from the standard Euler-Poincare reduction 
theorem, as provided in [14, 18]. We also have Sn t = S(gtr)t) = tu t + ^t^t- 
From ut = gtrjt, we get Sut = gtSrjt + Ct^t and from wt = gtSrjt we also have 
zb t = u t w t + g t rjt- This yields Su t = zu t + £ t u t - u t zu t . 

We also compute the boundary conditions for £ and zu. At t = 0, we have 
go = id and no = goVo = est which implies £o = and wq = 0. At i = 1, the 
relation girji = est yields £ini + cji = 0. 
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Now, the first variation is 
' '8£ • , . \ ( 8£ 



I ( I & - ad ^) + h + 6n <) + ( J£ I + ^ - n H ) dt = °- 

In the integration by parts to eliminate & and wt, the boundary term is 
((5£/6u)i \ £1) + ((<^/Ji/)i| wi). Using the boundary condition, the last term 
can be re-written 

We therefore obtain the boundary equation 



81.,. 81.,. 



0. 



The evolution equation for £ is 



d 8£ ,„ <M <5£ <M 



and the one for is 



5 ^ <M 81 

dt 8v 8v 8n 



We therefore obtain the system 



(5) 



d 8£ * 8£ 8£ 8£ 
— — + ad u — + — on t + — o v t 



dt 8u 



Ut 8u Sn 



8v 



d_M _ 8£_ 

dt 8v 1 8v 8n 

8£ . . 81.,. 

(1) + (!).„, =0 



k h t = vt + u t n t 

The system (5) is equivalent to (4), since they characterize the same 
critical points. A direct evidence of this fact can be obtained by rewriting 
the first equation in (5) under the form: 

(8£ 8£ \ ,„ ( 8£ 8£ \ n 

We indeed have, for a solution of (5), 
d_(8£_ 8£_ \ d_8£_ / d_M\ 8£_ . 

dt\8ut 8v V dtSu \dtSu) nt 8v Ut 



d 8£ ( 8£ 



8£ 



8£ 



9i 5n V <5n 5i/ 
9_« 8£_ M 

dt 8u 5n 1 8v 1 
- -ad,,.— -ad (— on t ). 



( 8£\ 8£ 



<5u 



5^ 
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In the last equation, we have used the fact that, for any a £ g, 

(jr- o (un) — (u* ^-J on = (J^ a(un) — u{an)^j 



/6£_ 

'81 



it, a\n 



- -(-on\[u,a}) 

( A* f 5£ 



This equation combined with (5i/5u)\ + (S£/Sv)i ou\ = obviously implies 
the first equation in (4). 



5. Special cases 

5.1. Riemannian metric. A primary application of this framework can be 
based on the definition of a Riemannian metric on G x iV which is invariant 
for the action of G: (g, rj)h = (gh, /i _1 r/), the corresponding Lagrangian then 
taking the form 

l(u,n,v) = \\(u,v)\\l. 

The variational problem is now equivalent to the computation of geodesies 
for the canonical projection of this metric from G x N onto N. This con- 
struction has been introduced in [20]. The evolution equations have been 
derived and studied in [28] in the case l(u,n,v) = \u\g + |z^, for a given 
norm, |.| , on g and a pre-existing Riemannian structure on N. 

The interest of this construction is that this provides a Riemannian metric 
on iV which incorporates the group actions. Examples of this are given below 
for point sets and images. 



5.2. Semi-direct product. Assume that iV is a group and that for all 
g G G, the action of g on is a group homomorphism: For all n,h £ N, 
g{nh) = (gn){gn) (for example, N can be a vector space and the action of G 
can be linear). Consider the semi-direct product G(S)N with (g,n)(g,h) = 
(gg, (gh)n) and build on G®N a right-invariant metric constrained by its 
value || ||(id G ,idjv) & t t ne identity. Then optimizing the geodesic energy in 
G(§)N between (idc^o) and (<7i,?ti) with fixed no and n\ and free g\ yields 
a particular case of metamorphoses. 

Right invariance for the metric on G(§)N implies 

\\(UgX)\\( g ,n) = \\(Ug , {Un)n + {gn)C\\(gg i( gn)n) 
which, using (g,h) = (g^ 1 , <7 _1 n _1 ), yields, letting u = Ug^ 1 , 

\\(Ug,0\\(g,n) = II («, {un' l )n + n^CH (idcidjv) 

= \\{u,n- l {C-un)\\ {iAGMN) 
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since = u(n l n) = (un l )n + n l {un). So, the geodesic energy on G(§)N 
for a path of length 1 is 



I \\(n,n t \fit-u t nt)\\l dGMN) 

J 



and optimizing this with fixed no and n\ is equivalent to solving the meta- 
morphosis problem with 



(6) 



l(u, n, v) 



1 ' / )ll(id G) id JV )- 



This turns out to be a particular case of the previous example. The 
situation is even simpler when N is a vector space since this implies rr x v = v 
for all n and the Lagrangian does not depend on n, which gives a very simple 
form to systems (4) and (5). They become 



(7) 



— + — on 
Su bv 







d 51 



51 



h t = u t + u t n t 



and 



(8) 



dt 5u 



a <Tu + Tu 01/4 = 



dt Sv Ut bv 

*; (1) + s;< 1) »" 1 



n t = u t + u t n t 



Even when N is not a vector space, metamorphoses obtained from the 
semi-direct product formulation are specific among general metamorphoses, 
because they satisfy the conservation of momentum property which comes 
with every Lie group with a right invariant metric. This conservation equa- 
tion can be written 



' SI 81 . 

iT'ir ] =Ad 

\0U ou 



(gt,m)- 



where the adjoint representation is the one associated to the semi-direct 
product. This property (that we do not explicit in the general case) will be 
illustrated in some of the examples below. 
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5.3. Constrained Metamorphoses. Returning to the general formula- 
tion, it is sometimes useful to include constraints on m, the image of the 
metamorphosis, making the minimization in (2) subject to = (i € 

[0, 1]) for some function $ : N — > R 9 . Using Lagrange multipliers, this 
directly provides a new version of (4), yielding 

se 

5u 



(9) 



se 

+ — o n t = 



dt bv Ut bv bn 1 bn 
h t = v t + u t n t 



( 5 ^\ ■ \ 







6. Examples from Pattern Matching 

In the following examples and in the rest of the paper, G is a group 
of diffeomorphisms over some open subset Q, C M. d . We will assume that 
elements of G can be obtained as flows associated to ordinary differential 
equations of the form g t = ut°gt where ut is assumed to belong, at all times, 
to a Hilbert space g of vector fields on Q with the condition 

l 

I I 2 T 

||wt|| fl cre < oo. 

We will assume that elements of q are smooth enough, namely that V can 
be continuously embedded in the space of C p vector fields with vanishing p 
first derivatives on d£l and at infinity, for some p > 1. More details in this 
construction can be found in [27] (Appendix C). 

We will write the inner product in q under the form (u, v} y = (L g u \ t>) 
where L g is the duality operator from g to g*. Its inverse, a kernel operator, 
will be denoted K g . 

6.1. Landmarks and Peakons. The space N contains the objects that 
are subject to deformations. The simplest case probably corresponds to 
configurations of Q landmarks, for which N = Q®. So elements r},v G H are 
Q-tuples of points in fi, with tangent vectors being Q-tuples of d-dimensional 
vectors. 

The model that has been proposed in [20, 5] corresponds to the Lagrangian 
l(u,n,v) = \\u\\l + ±Y,\» {k) \ 2 - 

° k=l 

This Lagrangian is therefore independent of n (but does not correspond to 
a semi direct product). We have 51/ 5u = 2L g u, and (S£/Su) = (2/a 2 )(i'^ , . . . . 



L 
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Let n = (q^\ . . . , q^). From the definition ((81 /8v) o n\ w) = — (8l/8v\ wn) , 
we get (since wn = (w(q^), . . . ,w(q^))): 

81 ^ (fe) x 

on = -2^ — ®8 q{k) . 

k=l 

Here and later, we use the following notation: if / is a vector field on R d 
(considered as a vector density), and a a measure on M. d , the linear form 
/ <8> //, acting on vector fields, is defined by 

(10) (f® a \w)= f(xfw(x)dn. 

Our first equation for landmark metamorphosis is therefore 

k=l 

The second equation is (8£/8v) +u*(8£/8v) = 0, which in this case gives 

4 k) +Dut(q^) T ul k) =0, k = l,...,Q. 
Introducing = /a 2 , we can rewrite system (4) in the form: 

Q 

L s u t = ^2Pt k) ® 
k=i 

< pi k) +Du t (q i t k) ) T p < t k) =0, k = l,...,Q 

, q^=u t (q { ' ) )+v 2 p ( i\ k = l,...,Q 
Putting the evolution equations into this form is interesting because the lim- 
iting case, a 2 = 0, exactly corresponds to the peakon solution of the EPDiff 
equation [14], the dynamics of which having been recently described in [19]. 
It is important to see, however, that the solutions may have significantly 
distinct behavior when a > 0. Figures 1 and 2 illustrate this in the case 
of two landmarks in ID. The plots show the evolution of r = q2 — q± over 
time when K s = L" 1 is a Gaussian kernel. Figure 1 provides a comparison 
in the case of a head-on collision (p\ + P2 = 0). In the case a 2 = 0, the 
peakons approach each other infinitely closely in time without colliding. For 
positive a 2 , the peakons get close, slow down, then cross over and their dis- 
tance grows rapidly to infinity. The duration of the collision phase (q ~ 0) 
decreases when the sum of the absolute momenta (p\ — P2) increases. Note 
that the deformation g t never becomes singular during this process. The 
space first contracts when the peakons get closer, then expands after the 
crossover. 

In Figure 2, the case of one peakon overtaking another is shown. In 
the case a 2 = 0, we observe a well-known behavior: the peakons approach 
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Figure 1. Head-on collision between two peakons. The 
plots show the evolution of q = qi — q\ over time for several 
values of p = P2 — Pi, for a 2 = (left) and a 2 = 10~ 4 (right) 




1 20 30 40 50 60 1 20 30 40 BO 



Figure 2. One peakon catching on with another. The plots 
show the evolution of q = q<i — qi over time for several values 
of p = P2 — Pi, for a 2 = (left) and a 2 = 0.05 (right) 

very closely, then separate again without crossing over. For metamorphoses 
(a 2 > 0), the details of the behavior depend on the initial difference between 
the momenta. If it is small, then the evolution is similar to the case a 2 = 
0. When the initial momentum difference becomes larger, the peakon that 
started behind has enough energy to overpass the other one and the two 
peakons exchange position. In all cases, the deformations first experiences a 
contraction, then an expansion (relative to the position of the two peakons). 

6.2. Images. Now, consider the case when N is a space of smooth functions 
from SI to M, that we will call images, with the action (g,n) t— ► n o g . A 
simple case of metamorphoses [20, 27] can be obtained with the Lagrangian 

l{u,v) = \\u\\\ + \\\v\\ 2 L 2. 
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If w £ q and n is an image, wn = —Vn T w, so that (J£ on| w) = (|^| Vn T w) . 
Thus, since 5i/5v = 2u/a 2 , the first equation is 



LgU t = 

Since u* (5£/5v) is defined by 



; u t Vn t . 



(51 N 



Vw T ?t 



--(^IV^u) 
^2 (div(i/u) | , 



we obtain the second equation 



^ + ^2div(^ t u t ) = 0. 

As in the landmark case, denote z = v/a 2 and rewrite the evolution equa- 
tions in the form 



(12) 



L s u t = -z t Vn t 
z t + div(z t u t ) = 



fit + Vnfut = a 2 z t 



Existence and uniqueness of solutions for this system have been proved in 
[27]. From a visual point of view, image metamorphoses are similar to 
what is usually called "morphing" in computer graphics. The evolution of 
the image over time, t t— > nt, is a combination of deformations and image 
intensity variation. Algorithms and experimental results for the solution of 
the boundary value problem (minimize the Lagrangian between two images) 
can be found in [20, 9]. Some examples of minimizing geodesies are also 
provided in Figure 3. 



In ID, letting m = L s u = (1 — d 2 )u, the time evolving form of this 
system (as provided by (5), or by direct computation from (16)) becomes, 
with p = az: 

(13) dtm + ud x m + 2md x u = —pd x p with dtp + d x (pu) = 

This relates, with the important difference of a minus sign in front of pd x p 
in the first equation, to the two-component Camassa-Holm system studied 
in [6, 8, 17]. The system (13) in our case is equivalent to the compatibility 
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W W W * « i* 



Figure 3. Minimizing metamorphosis between images. 
The optimal trajectories for m are computed between the 
first and last images in each row. The remaining images 
show ni at intermediate points in time. 



for dX/dt = of 

(14) 

(15) 



+ 



- + m\ + p 2 \ 2 )ip = 
+ ujd x ip + -ipd x u 



Image matching can also be seen under the semi-direct product point of 
view, since the action is linear and the Lagrangian takes the form (6) with 
n^v = v. This implies that the momentum, which is, in this case, the pair 
(L g u, z), is conserved in a fixed frame. Working out the conservation equa- 
tion Ad* 5 n j(L g ii, z) = est in this case yields the equations L s ut+ztVn t = est 

and zt = det(Dg^ 1 )zQ o , This last condition is the integrated form of 
the second equation in (12), while the first equation in (12) implies that in 
fact L s ut + ztVnt = 0, which is the horizontality condition in the quotient 
space G®N/G. 

6.3. Densities. We here let N be a space of smooth functions n : — ► R 
with the action (g,n) i— > | det D(g~ )\n o g^ 1 , i.e., n deforms as a density. 
We consider the same Lagrangian as with images, 



£(u,u) 



1, 



1 



For w € 5 and n £ H, we have wn 
implies 

— onto = 
ou 1 



yielding the first equation 



-Vn T w — ndiv(u>) = — div(mu). This 



-— 2 {y \ div(nu>)) 



a- 



T riVv. 



a* 
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o 
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o 
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Figure 4. Minimizing metamorphoses between densities 
with equal total mass. The optimal trajectories for nt are 
computed between the first and last densities in each row. 
The remaining images show m at intermediate points in time. 

Similarly, we get u~kv = Vv T u and the equation 

v + Vv T v = 0. 

This yields the system, where we have, as before, introduced z = vja 2 : 
(16) 



z + Vz T u = 



a 2 z 



n + div(nu) 

We are here also in the semi-direct product case, the equations for the 
conservation of momentum being L g u + nVz = est and z = z^og^ 1 . Like for 
images, the constant in the first conservation equation vanishes for horizontal 
geodesies in G©N/G. Optimal metamorphoses with densities are illustrated 
in Figure 4. 



6.4. Plane Curves. We here consider matching unit-length curves 7 de- 
fined on the unit circle S , represented, as in [32, 21, 33] with their nor- 
malized tangent 9 1— > 70 with |7e| = l/2n. The set N is therefore a set of 
functions n : S 1 — > S' 1 (l/27r) where 5' 1 (r) is the sphere with radius r in M 2 . 
We let G be the group of diffeomorphisms of 5 1 and consider the reduced 
Lagrangian 

£(u,v) = I u 2 e d6 + \ [ \is\ 2 d6. 
Js 1 & Jsi 

We want to solve the metamorphosis problem while ensuring that curves are 
closed, which translates into: 

-1 

n t d9 = 0. 



To explicit (9), we need a local chart to compute the partial derivatives 
d/dv and d/dn (N is not a vector space here). Consider the representation 
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Figure 5. Optimal metamorphoses between plane 
curves. The optimal trajectories for rit are computed between 
the first and last curves in each row. The remaining images 
show n t at intermediate points in time. 



n = h a and v = a 2 ph^ with h a = (cos a, sin a) and 
We then get the equations, with A ( 6 M 2 , 

d 2 ut dott 

~m 2 ~ + pt ^0 =o 



(- 



sm a, cos a 



(17) 



dp t d 

+ de {utpt) 



<y 2 Pt 



Uf 



dot 
dd 



^ Js 1 



uth^dO = 



Interestingly, these equations can be notably simplified in the case a 2 = 1, 
which has been considered in studied in [32, 25, 26, 33]. In this case, the 
change of variables z 2 = gtntogt, where both sides are interpreted as complex 
numbers reduces (17) to the geodesic equations on a Grassmann manifold, 
on which explicit computations can be made [33]. The case a 2 = 4 is also 
interesting, and has been discussed in [15]. Figure 5 provides the result of 
curve metamorphosis using a 2 = 1. 



7. Measure Metamorphoses 

We now focus on extending the example in section 6.3 to include also 
singular measures. The L 2 norm that we have used between densities is 
therefore no longer available. We will rely on a construction that was intro- 
duced in [10]. 

We let H be a reproducing kernel Hilbert space (RKHS) of functions 
over M. d and N = H*; H being an RKHS is equivalent to the fact that 
for all x € M d , there exists a function K x G H such that, for all / G H, 
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(K x , f) H = f(x). The kernel Kjj{x,y) := K x (y) satisfies the equation 
(K x , K y ) = Kh(x, y). It also provides an isometry between N and H via 
the relation 77 ^ Kur] with (Kht] , f} H = (r)\ /). This implies in particular 
that the dual inner product on N is given by 

i*!' = H K HV)- 

Letting G be like in the previous section, we want to study metamorphoses 
in G x N. We define the action of G on N by (with g G G,r) G N and / G H) 

(av\ /) = {v\f°9)- 

This obviously generalizes the action on densities discussed in the previous 
section. For 77 G N and w G 0, we have, for all / G H, (wrj\ /) = Vf T w). 

Since the action is linear, we can use the semi-direct product model and 
the Lagrangian 

£(u, v) 



1 



To explicitly compute again (4) in this context, we need to compute 
(dljdv) o n. We let / = bijbv = (1/<j 2 )Khv- By definition, we have, 
for all w G 



(/ o n\ w) 



"(/I 



- In 
(V/ 



wnj 
n\ w) 



Since (/| urj) 



We therefore obtain our first equation (5t/5u) = V/ £2 

Vf T u), the second equation, ft + ut*ft = is the advection: f t + 
V T /t«t = 0. This yields the system (with Lh = K^ 1 ) 

L s u t = V/t <8> n t 



(18) 



A + v/ t r ^ = 



k n t - u t n t = a 2 L H f t 
From the second equation, we get ft = fo° 9t l ■ Using gr\ 



<J 2 L H f 



and grj = n, we get 

n t 



g t n + a 2 g t I g s 1 K H 1 (f o s- s ^ds. 

JO 



We therefore obtain integrated equations for measure metamorphoses 

L s u t = V/t <8> n t 
(19) I - t 

\ Jo 

with g t = u t og t . 
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We now pass to the theoretical study of the existence of solutions for the 
initial value problem (IVP) and boundary value problem (BVP) for measure 
metamorphosis (with uniqueness in the IVP case). The next two sections 
are notably more technical than the rest of this paper. They are well isolated 
from it, however, and it is possible, if desired, to skip directly to section 10. 

8. Existence of solutions for the measure metamorphosis IVP 

8.1. Hypotheses on the Hilbert spaces. For the existence proofs to 
proceed, we need some conditions on the Hilbert spaces g and H. They are 
essentially adapted to H being equivalent to Hq (the completion, in H q , 
of C°° functions with compact support), in which case H~ = H^' 1 and 
H + = Hq +1 (q being large enough to ensure that H q is embedded in C°). 

In the following, est represent a generic constant, and C a generic contin- 
uous function of its parameters. We assume the existence of two spaces H + 
and H~ and the following properties, valid for some q > 1. 
(HI) If H = H~, H or H+, we have: if / G H and g G C q (Q), then fog£H 
and 

\\fog\\ 6 <C(\\g\\ qt0O )\\f\\ 6 . 
(H2) For / G H + , and g,g two C q diffeomorphisms 

\\f°g-f° 9\\h < cst||/|| H +C(max(||fli|| g)00 , ||ff|| g ,oo))||ff - <7||?,oo- 

(H3) For / S H, define the operator Qj on q by Qjw = X7f T w. Then for 
all / G H and g G C q (Q), Qf maps X q {Q) to H~ with 

\\Qfw\\n- < cst ll/ll 

for all w G X q (Q). (Here, X q (VL) denotes the set of C q vector fields on Q 

with the supremum norm of all derivatives of order less than q.) 

(H4) If / G H+, then K~ l f G {H~)* and for all z G H~, {K H l f\z) < 

cst\\f\\ H +\\z\\ H -. 

(H5) g is continuously included in Xq (fi) for p > q + 1, where Xq is the 
completion of compactly supported vector fields in X P (Q). 
Denote N + = (H~)*. We then have 

Theorem 1. Under hypotheses (HI) to (H3), for all T > 0, there exists a 
unique solution to system (19) over [0, T] with initial conditions uq G N + 
and fo G H + . 

Proof. We prove existence for small enough T with a fixed point argument. 
Consider the Hilbert space L 2 ([0,T],g), with norm 
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For u G L 2 ([0,T],g), define := u' given by 
( u' t = K 8 (Vf t ®n t ) 

(20) | ft _ 

n t = g t n + a 2 g t J g s 1 L H (f o g s 1 )ds 

with g t = u t o g t . 

First note that the hypotheses imply that ^ is well defined and takes 
values in L 2 ([0, T],g). Indeed, by definition, 

(u' t ,w) g = {Vf t ®n t \w) 

= (n t \Vf?w) 

< cst||nt||jv||/t||ff+|Mkoo 

< cst||nt||Ar||/t||H+lk|| 8 

so that a sufficient condition for V/t<8mt G 0* is n t G N and ft G Since 
ft = foogi 1 , we have \\f t \\ H + < est || / 1| j^-h Cdl^ 1 IU,oo)- From [27], we have 
H^r 1 IU,oo = 0{\\u\\ 2 ,t), yielding 

(21) ||/t||H+<CSt||/o|| H+ C(||«|| 2) T). 
For z G H, we have 

(n t \z) = (n \ z o g t ) + a 2 [ (KJj 1 f s \z s o g t o g'^ds 

Jo 

||no|kC(|| 5t || (? , 00 )+a 2 ^ HMI^Cdl^o^- 1 !!^)^ \\z\\ H 
so that 

(22) IKHjv < (UnolU + ^H/oIIh) C(||«|| 2 ,t). 
This implies that 

Kilo < ll/o||if+ (llnolk + ^ll/ollif) C(||u|| 2 , T ) 

and 

Hu'lla.T < cstVr||/o|| H + (IKIk + <T 2 T\\fo\\ H ) C(\\u\\ 2 ,t). 
In particular, this implies that, for any M > 0, there exists a Tq{M) (only 
depending on ||no||jv+ and ||/o||h+) such that, for T < To, ||«|| 2 ,t < ^ 
implies ||u'|| 2 ,t < From now on, we assume that T < Tq(M) with 

M = l. 

Note that a similar computation also shows that n t G iV + : for z G 
we have 

(n t \z) = (n \ z o g t ) + a 2 (K^ 1 f s \z s o g t o gj^ds 

Jo 

< (\\no\\N + C(\\g t \\ q ,oo)+cT 2 ^ \\fs\\H+C(\\g t og; l \\ qt00 )ds^ \\z\\ H - 



< 
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so that 

(23) K|U + < (IKHJV+ +<7 2 *||/o||h+)C(||u|| 2) t). 

We now ensure that ^ is contractive. Take u, u with max(||u||2,T) IN|2,t) < 
1. We want to show that T can be chosen so that \\u' — u'\\2,t < P\\ u ~ ^l|2,r 
with p < 1, where u' = ^f(u) and u' = ^(u). We have 

IK -"tils = l|V/t<8>n t - V/ t (8)n t || r 

< ||V(/ t - ft) ® n t || r + \\Vf t ® (nt - nt)\\g> 

< cstllntlljv+H/t - f t \\H + est || nt - nt||jv||/t||if+ 

(Here, we have used the fact that for n G N + and z £ H, we have both 
(n| z) < ||n||jv||z||ff and (n| z) < ||n||jv+ IMIif- ■) 

Upper bounds for ||^||jv and ||nt||jv-+ are provided by equations 

(21), (22) and (23). 

Moreover, from (H2), we have 

||/t-/t||jf < CCmaxdl^H^oo.ll^ll^ooJII/ollH+ll^" 1 -g^Koo 
and, using [27], we have ||<7 f _1 — g^ 1 \\ q +i t oo < cst||u — tt||2,T so that 

(24) KlU+ll/t-ZtllH <™t(\\n \\ N+ +ah\\f \\ H+ )\\f \\ H+ \\u-u\\ 2 , T 
For z G H, we can write 

(n t -n t \z) 



(i) 
(ii) 
(iii) 

For (i), we have used 

(n 1 zog t - zog t ) 

< \\no\\N+\\ z ° 9t ~ zog t \\ H - < cst||n ||iv+||2:||fl-||gt - <7t|k<x>- 

For (ii), we have used the same argument combined with the fact that, 
since f a € H+, K' 1 ^ G N+ with \\K^-f 8 \\ N+ < cst||/ a || H +. For (iii), the 
computation uses the fact that for z G H, 

(K H Hf s -f s )\z) < II/.-AHhIIzIIh. 

We therefore obtain the inequality 

(25) \\n t - h t \\ H < cst(||n ||j V + + ||/o||h+)||« - "Ikr- 



= (n \ zo g t - zog t ) 

+ " 2 / { K H l fs\ zo 9t° 9s l - z°9t° gl^ds 
Jo 

W I (Ka\fs - f»)\z°9t°97 X )d* 
Jo 

< cst\\n \\ N +\\z\\ H \\gt - gt\\q,oo 

+ est a 2 t\\f \\ H+ \\z\\ H max\\g t o g' 1 - g t og^ l \\ q ,oo 

+ CSt cr 2 t||/ ||//+||z||^ max||5- s - g s \\q,oo\\9t° SiJ^Woo- 
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Collecting the previous estimates, we have 

IK - KWa ^ F (\\ n o\\N+, ||/o||h+)||« - «||2,t, 
where F is a polynomial. This implies 

\\u' - u'\\ 2 ,t < VTF(\\n \\ N +,\\fo\\ H +)\\u - uh,T, 

so that \I> is contractive for small enough T. 

The extension from small times to all times can be done like in [27], and we 
only sketch the details. According to the small time computation, the length 
over which the solution exists is at least the inverse of a polynomial function 
of ||no||jv+ an d ||/o||ff+- One will therefore be able to extend this equation 
beyond T, unless either ||nt||jv+ or tends to infinity when t tends to 

T. From (21) and (22), this can happen only if ||u||2 1 tends to infinity when 
t tends to T. But this is impossible, because (19) is a geodesic equation on 
G®N which implies that the value of h t := \\ut ||g + (1/<t 2 )||^||^ is constant 
over time. This implies in particular that 1 1 1 1 2,* — ^0 and therefore cannot 
tend to infinity in finite time. □ 

9. Existence of solutions for the measure metamorphosis BVP 

Our goal in this section is to prove that, under some conditions on H, 
the boundary value problem for measure metamorphoses (BVP) always has 
solutions. This problem requires to minimize, with fixed uq and n\, 

E(u,n) := I \\u t \\ 2 dt + ^ f \\h t - u t n t \\ 2 N dt. 
Jo ' Jo 

Letting f t = Kjj(h t — u t n t ), the problem is equivalent to minimize 

E(u,f) := f 1 \\u t \\ldt + a 2 C \\f t f H dt 
Jo ' Jo 

with boundary 

n\ = gi n + a 2 gi I gj 1 ^^ 1 f s )ds. 
Jo 

We have 

Theorem 2. Assume that (HI) and (H2) hold and that H + is dense in H. 
Then, for any given no,ni G H, there exists a minimizer for the BVP. 

Proof. Consider a minimizing sequence (u^ k \ f^). Using a time change 

if necessary, we can ensure that := ||?4^||g + o^H/t^H 2 is indepen- 
dent of time and therefore bounded. Moreover, we can extract a subse- 
quence (still denoted (u^ k \ /^)) which weakly converges to some (u, f) in 
L 2 ([0, 1], g x H), equipped with the norm ||(it', _f')|| = E(u' , /'). This implies 
that E(u, f ) < liminf E(u^ k \ f^) so that the only thing that needs to be 
showed is that the boundary condition is still satisfied, namely, letting 

n' 1 =gin + a 2 g 1 g^ l (K^ l f s )ds 
Jo 



20 DARRYL D HOLM, ALAIN TROUVE, AND LAURENT YOUNES 

with g t = ut o gt, we have n[ = n\. For this, it suffices to show that, for z 
in a dense subset of H , we have (n' : | z) = (ni | z) . Since, for all k, we have 

Jo 

it suffices to show that, for z in a dense subset of if, we have 

(26) (n 1 z o -> (n 1 z o 51) 
and 

(27) r (^Vi fc) | z o fl <*> o CgW)" 1 )* - f {K-'f s \ zo 9l o g^)ds 
Jo Jo 

Because of the weak convergence of to u, the flows g^ converge to 
g for the (p — 1, oo)-norm, uniformly in time [27]. Taking z £ we have 

(n , zojf - zo 9l )< 

C(sup(||s' i (fc) ||p,oo, ||5t||p, oo ))||n ||fl-*||2:||jj-+ sup - fft||p,oo 
t t 

which tends to so that (26) is true for all z £ H + . Splitting the terms in 

(27) , it suffices to show that for all s 6 [0, 1], 

(28) /V*Vi fe) |*°#M#r>- fW I "floj-VO 

JO JO 

and 

(29) f {K H l f^\z o 5l o {g s y l )ds - f 1 {K H 1 f s \ zo 9l o gj^ds - 0. 

JO JO 

For (28), and because we have ensured that ||/i*^||.ff and ||w^|| are uni- 
formly bounded, it suffices to show that each term in the integral tends to 
and then use the dominated convergence theorem. The left-hand term of 
(28) is bounded in absolute value by 

cst\\fW\\ H \\z\\ H+ \\g[ k} o ( 5 ( fe ))-i _ gi gJ% >00 

which tends to 0, so that (28) holds. For (29) we only need the fact that 
z 9i 57 1 belongs to H and the weak convergence of to / to conclude. 
This shows that n[ = n\ and concludes the proof of the theorem. □ 

9.1. Remark. Equations (19) have been obtained from general formulae 
that were derived under the assumption that G is a Lie group, (which is 
not the case here). It is important to rigorously recompute the Euler equa- 
tion to reconnect the IVP and the BVP. The variation with respect to u is 
straightforward and provides the first equation in (19). 
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We now discuss the minimization in n with fixed u. Letting t] = g~ l n the 
problem with fixed u (and therefore fixed g) is equivalent to the minimization 
of 

F(rf) = / II^Hat^ 
J o 

with fixed boundary conditions t]q = uq and r\\ = g^ni. 
Assume the following hypothesis: 

(Hlb) z — > z o g is weakly continuous for z £ H (where H* = N) 

For g £ G, introduce the operator K 9 H defined by K 9 H r] = Kij{grj) o g. 
With this notation, 




Let L 9 H = (K 9 H ) 1 , i.e., L^/ = p 1 .L ff (/ o g 1 ), and denote 

(30) L 9 H = f L 9 *dt. 

Jo 

We have for any / £ H, 

(31) (iy|/)= f\l%f\f)dt= f'wfog^fHdt^cwfwl 

Jo Jo 

where the last inequality comes from (HI) with C = [sup t Cdl^ 1 1| ?;00 )] _1 . 
Thus, L 9 H f = if and only if / = and L 9 H {H) is dense in H*. If we 
prove that L 9 H (H) is closed, we will get that L 9 H (H) : H — > H* is invertible. 
Let (f n )n>o be a sequence in iT and rj £ H* such that L 9 H (f n ) — > 77 in 
Then we get from (31) that / n is bounded in H and we can assume that it 
admits a weak limit in H . Thus (£#/oo| /) = / (/oo 5t~\ f°9t 1 )Hdt = 
lim /^(/n ° / ° 9^ 1 }ndt where the last equality comes from (Hlb) and 
the dominated convergence theorem. This yields {L 9 H foo\ f) = /) f° r 
any f £ H so that 77 = L 9 H fo2 and L 9 H (H) is closed. 
Let /o such that ^/o = Vi ~ Vo an d define 

(32) r?t = r?o+ f L 9 £fods 

Jo 

for any i G [0, 1]. If we can prove that for any (rf t ) £ L 2 ([0, 1],H*) with 771 
and 770 fixed 

F(rj — 77) = i^r/) + est 
we deduce immediately the result. However, we have 

F(ri -fj) = F( V ) + F(r?) - 2 (r) t | 

Jo 

and /J (r) t | K 9 ^fj t )dt = J* (r)t\ fo)dt = (771 - 770 1 /o) so that the result is 
proved. 
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We therefore have proved that minimizing E with respect to u, n is the 
same as minimizing 

E(u)= / \\u t \\ 2 dt + —(n 1 -r] \(L g H ) 1 (m~Vo)) 
Jo ° 

with respect to u. We also have the expression (32) for the optimal n which 
is consistent with the second equation in (19). 

10. A COMPUTATIONAL ANSATZ FOR POINT MEASURE MATCHING 

The previous theorems provide a rigorous foundation for the considered 
measure matching approach. However, an important issue needs to be ad- 
dressed. The space N, which has been introduced in order to take advantage 
of its Hilbert structure, is a very big space that contains distributions that 
are more singular than measures. Now, when matching two measures no 
and n±, the question naturally arises of whether the optimal evolution, i.e., 
the measure nt, can turn up being more singular than measures, since the 
existence theorem only ensures that it belongs to N. 

The second equation in (19) indicates that this should not be the case, 
since it says that 



n t = gtm + 9t [ l h( l h) 1 (.9i ln i ~ no 
Jo 



)ds 



where L g H and L 9 H are defined above and in equation (30). Thinking of 
Lh as a differential operator, the number of derivatives computed by L 9 ^ is 
"canceled" by the pre-application of (L 9 H )~ l so that nt should not be more 
singular than (g^ 1 n± — no) which is a measure. It is therefore reasonable to 
conjecture that when no and n\ are weighted sums of Dirac measures (which 
is a case of practical interest), then n t is a measure which has an absolutely 
continuous part, and a singular part which is also a sum of Diracs. More 
precisely, if 

q r 

k=i k k=i k 

a reasonable ansatz for nt is 

r r 

n t = Y^ a k(t)5 Xk (t) + Pk(t)S yk (t) + f(t, .)dx. 

k=l k=l 

Assuming this, the Lagrangian £(vt,n t , v t ) can be considered as a function 
of («fc(t)), (/?fc(i)), (xfc(t)), (yk(t)), f(t, •) and their time derivatives, with an 
explicit expression in terms of the kernel Kh and its space derivatives that 
we do not provide here, since it is quite lengthy. Minimization can then 
be done with standard methods, with boundary conditions a&(0) = ot^\ 
P k (0) = 0, a k (l) = 0, /3 fc (l) = P { k 1] , x fc (0) = xf- I/fc(l) = Vk ] and /(0, •) = 

/(i, -) = o. 
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11. More Metamorphosis 



Without getting into the level of rigor and detail developed with mea- 
sure metamorphosis, we now review additional situations in which meta- 
morphoses can be used. At the difference of the examples considered in 
section 6, the following models have not been solved numerically yet, nor 
has a theoretical analysis been developed although we expect that the pre- 
vious proofs of existence of solutions can be modified to work in these cases 
also. 

11.1. Singular Image Metamorphosis. In section 7, we have extended 
density metamorphosis to a context that includes singular measures. A 
similar analysis can be made to extend image metamorphosis to generalized 
functions. Let H be a space of smooth scalar functions and N = H* as 
before. The extension to N of the action of diffeomorphisms on images is 
(g, n) ' ^ gn with 

(gn\ f) = (n| det(Dg)f o g) , 

the infinitesimal action being (un\ f) = (n| div(/w)). Take as before the 
simplest reduced Lagrangian 

1 



l{u,v) = \\u\\ 2 s + -^\\u\\ 2 N . 



We have, letting / = (1/o 2 )Khv 

,se 

5u 



,51 , , 
(j7, on \ u ) 



At 51 \ \ 

and [u -k — u3) 

du 1 ' 



(/Vn| u) 
(uj\ div(/u)) 



with the notation (Vn| w) = — (n\ divw). We therefore obtain the general- 
ized version of (12): 



(33) 



L s u t = -f t Vn t 
ft + dw(f t u t ) = 



k h t - u t n t = a L H f t 
Like for measures, this leads to integrated equations 
L g u t = -/tVn t 



(34) 



rt 

n t = gtn + cr 2 gt / g s 1 L H (/ o g s 1 det(g s 1 )) ds 
K JO 



with g t = u t og t . 

In the ID case, with p = af and m = L Q u = (1 — d 2 )u, we get a new 
version of (13): 

(35) dtvn + ud x m + 2md x u = —pd x Lnp with dtp + d x {pu) = 
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11.2. Smooth Image Metamorphosis. We can go in the opposite direc- 
tion and consider norms that will apply to smooth images in the metamor- 
phosis formulation. Namely, keeping the notation of the previous section, 
we can define 

1{U,V) = |H|g + — \\V\\ 2 H . 



Since the situation is completely symmetrical, we can immediately write the 
new system, with / = (1/<7 2 )Lhv, 



(36) 



L s u t = -ftVnt 
ft + dw(f t u t ) = 
h t - u t n t = a 2 K H f t 



The second equation must be understood in a generalized sense, / being 
advected by the flow as a measure. 

An interesting feature of this system is that it admits singular solutions 

for the pair (L s u, /). Indeed, assume that /o = Ylk=i (*) is a sum °f 

x o 

weighted point masses. The second equation in (36) implies that 

k=i 



with x[ k ^ = u t (xf'). We also have 



Q 



L S U t = Yl ^ S : 



.(fc) 



k=l 



with a^' = w^Vn(x t ) and n evolves according to the last equation in the 
system. 

Exploiting this very simple structure can lead to interesting new methods 
for the analysis of smooth (or smoothed) images and will be considered in 
future work. The ID evolution equation associated to this context is (letting, 
again, p = af) 



(37) 



dtvn + ud x m + 2md x u = —pd x Knp with dtp + d x (pu) = 



11.3. Smooth Densities. Obviously, the same reduced Lagrangian can be 
used with densities. Using (V/| uj = — (f\ div(it)) as a definition of V/ for 
generalized functions, we get the system: 



(38) 



L s u t = -Vntft - ntVft 
ft + nfVf t = 
k ht + div(u t n t ) = a 2 K H f t 
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Here also, singular solutions in (L g u, f) exist, although (as seen from the 
first equation), L s u now is one derivative less regular than /. 

12. Conclusion 

We have provided a general framework for the pattern matching theory 
of metamorphoses. The provided equations are quite versatile, and adapted 
to any context in which a Lie group acts on a manifold. 

In the particular case of diffeomorphisms acting on generalized functions, 
we have obtained a new set of equations, and showed that they had solutions, 
and so does the initial variational problem. Our equations seem to indicate 
that, when matching measures, metamorphoses do not generate additional 
singularities, but that they may introduce smooth components that did not 
appear in the initial problem. 

Open for future work is the interesting problem of building an efficient nu- 
merical implementation of measure metamorphoses, and their use in specific 
pattern matching applications. 
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